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We investigate the structure and collective modes of a planar surface of a trapped Bose-Einstein 
condensed gas at zero temperature. In the long-wavelength limit we find a mode similar to the 
gravity wave on the surface of a fiuid with the frequency uj and the wavenumber q related by 
uj^ = Fq/m. Here F is the force due to the confining potential at the surface and m is the particle 

. . . mass. At shorter wavelengths we use a variational approach and find corrections to tj^ of order 

O^ ' g^lng. We demonstrate the usefulness of the concept of an effective surface tension for describing 

0^ ' both static and dynamic properties of condensed atomic clouds. 
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I. INTRODUCTION 



Properties of trapped clouds of Bose-Einstein condensed atoms have been investigated intensively both experimen- 
tally and theoretically over the past few years, since the experimental realization of Bose-Einstein condensation in 
^ ] dilute atomic gases [Q. At zero temperature, the behavior of the order parameter, ijj, is determined by the time- 
dependent Gross-Pitaevskii equation, a Schrodinger equation with the nonlinear term proportional to IV'PV' added: 



^"^ ■ 2 

CN' _^V2^(r,t) + T/(r)V'(r,t) + ;7o|V(r,t)P^(r,t)=*?i^^^. (1) 

c_^ ' zm at 

Qv Here Uq is the effective two-particle interaction, which may be expressed in terms of the scattering length a according 

■^ to Uq = Airah jm, where m is the mass of an atom. The confining potential is denoted by Vij). The static structure 

C^ of the cloud is determined by the time-independent Gross-Pitaevskii equation, which reads 
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^ . - ^W-W + ^(r)V^(r) + t/o|V(r)rV(r) = ^i^iy), (2) 

c3 : . ™ 

Q ^ where \i is the chemical potential. For clouds containing a sufficiently large number of atoms with repulsive interactions, 
O ' many properties may be calculated to a good approximation using the Thomas-Fermi approximation [g| , in which one 
neglects the kinetic energy term in the Gross-Pitaevskii equation. One then finds a density profile 

X' |^(r)|2 = (A.-y(r))/C/o, /i>^(r), IV-WP^O, M<nr)- (3) 

H ; 

Cu The Thomas-Fermi approach may also be used for calculating collective modes, and this has been done by Stringari 

A number of phenomena are associated with the surface region of the cloud. One is that the contribution of the 
kinetic energy term to the total energy of the cloud comes mainly from the surface region. Also Stringari H has 
identified surface modes of oscillation in the Thomas-Fermi approximation. To understand such phenomena in detail, 
it is useful to consider the properties of a planar surface, and to approximate the trapping potential by a linear 
function of the coordinates. The potential is then given by 

V{x) = Fx, (4) 

where the coordinate x measures distances in the direction of VF. The presence of the V^ term in the Gross-Pitaevskii 
equation leads to a rounding-off of the density profile in the surface region over a distance of order 



n^ ^ '^' 



2mF 



(5) 



and within this ap£roach the static structure of the surface region has been studied, and contributions to the total 
energy calculated 



The purpose of this paper is to consider properties of plane surfaces of Bose-Einstein condensed clouds in a linear 
potential. This is of interest both for giving analytical expressions for a number of properties of large clouds, as 
well as for giving physical insight. First, we consider the static properties of the surface region, and show how the 
kinetic energy term in the Gross-Pitaevskii equation gives rise to an effective surface tension. We also consider the 
density profile in the vicinity of the surface, and show that the Gross-Pitaevskii equation and the Thomas-Fermi 
approximation to it lead to identical results for the column density of atoms down to a point well within the cloud. 
In the second part of the paper we consider surface collective modes. We begin by considering the Thomas-Fermi 
approach, which yields a mode with the dispersion relation w^ = {F/m)q, where w is the angular frequency and q the 
magnitude of the wave vector for the mode. This is exactly the same form as for gravity waves on the surface of a 
fluid, with the role of gravity being played by the trapping potential. At higher wave numbers there are contributions 
to the mode frequency arising from the V^V term in the Gross-Pitaevskii equation, and we evaluate the effect of these 
by employing a trial wave function that allows for the rounding-off of the order parameter profile in the surface region. 
The extra contributions to uP' are proportional to q^ ln(l/g(5), and they will be shown to have a ready interpretation 
in terms of the effective surface tension introduced in the discussion of static properties. 

II. STATIC PROPERTIES 

In this section we introduce the concept of an effective surface tension, and consider the number of particles that 
are associated with the rounding-off of the density profile in the surface region. We shall take the origin of the 
x-coordinatc to be at the point where the potential is equal to the chemical potential, which is where the density 
vanishes in the Thomas-Fermi approximation. 

Let us now calculate the contribution to the energy coming from the kinetic energy term in the Hamiltonian. As 
was discussed in Ref. B, the total kinetic energy may be written in terms of either ^_^ ("0') dx or j_^ ^pip"dx, where 
the prime denotes differentiation with respect to x: the kinetic energy density is not a well defined quantity. For 
calculating the total energy either expression may be used, but it must be used consistently. In the case of the surface 
problem under investigation, we wish to be able to associate a contribution to the kinetic energy with the surface, and 
there is no unambiguous way of doing this, since one could use either of the energy expressions integrated to some 
point well inside the surface (rather than — oo as in the case of the total kinetic energy). However, these differ by a 
constant. For definiteness, we shall employ the symmetrized expression and define the kinetic energy per unit area of 
the surface as 



^^=I^L(^')'^^' 



(6) 
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where the point x = —L is chosen to lie well inside the surface. From Ref. B this may be seen to be given by 

n^ F ^f L \ 
8^C^ 



where the coefficient 0.240 was found by numerical integration. If we had chosen the other expression for the kinetic 
energy density, the result would differ by h^F/imUo. The coefficient of InL is independent of the choice of kinetic 
energy density, but terms independent of L are not. This difference, however, is unimportant in application to surface 
modes, where the logarithmic term is the one of most interest. If this expression were independent of L, the kinetic 
energy would have precisely the same form as a surface tension, since the energy would be proportional to the area 
of the surface. Because the kinetic energy density falls off only slowly away from the surface, the total kinetic energy 
increases logarithmically with L, but this dependence is sufficiently mild that it is still useful to think of the kinetic 
energy as being physically analogous to a surface tension. 

Next we consider the density profile, and compare the Thomas-Fermi result with the exact one. For \x\ 3> 6 the 
Thomas-Fermi wave function is a good approximation to the exact one and thus tp" is negative. Consequently, as 
one can see from the Gross-Pitaevskii equation, the density is depressed below the Thomas-Fermi value. Outside the 
cloud, the density is increased, due to the quantum-mechanical tail of the wave function. We shall now demonstrate 
the somewhat surprising result that the number of particles associated with the depression of densities inside the 
cloud is exactly equal to the number associated with regions of increased density. 

To prove this result it is simplest to consider the momentum density, g, defined by 

. = A(W;*-^-W,^) (8) 

2i dxi dxi 



whose time dependence may be determined from the Gross-Pitaevskii equation (|i|). By differentiating (|8|) and using 
(|l|) we find that the momentum density satisfies the conservation condition 

dg, dT,j _ dV 
dt dxj dxi ' 

where n(= IV'P) is the density. Here the stress tensor Tij is given by 

T,,^T^j+P5^j, (10) 

where 

^■' 2m \dxi dxj dxi dxj 2 dxidxj J 
is the free-particle stress tensor, while P is the pressure arising from the interaction, 

P = \n^Uo. (12) 

For the one-dimensional problem we are considering, the condition for hydrostatic equilibrium is simply 



dx 

from which it follows that 



= -Fn, (13) 



ndx, (14) 

since the stress tensor vanishes for large positive x. Deep within the cloud the order parameter varies as (— x)^' ■^, and 
therefore T^^{—L) ~ 1/L, which vanishes for L ^ oo. Consequently 

' 1 

ndx^-7?{-L)Uo + 0{l/L). (15) 

The first term is the result one obtains in the Thomas-Fermi approximation, where n — —xF/Uq for cc < 0, and thus 
one sees that the total column density outside a point well within the surface is the same as that in the Thomas-Fermi 
approximation, apart from corrections of order 1/L. Physically the result follows from the condition that the stress 
tensor at any point must be balanced by the force on all the material at larger values of x, a result familiar in the 
context of equilibrium of fluids and gases in the presence of gravitational fields. 

III. SURFACE MODES OF OSCILLATION AT LONG WAVELENGTHS 

We turn now to time-dependent situations and consider oscillations. For this purpose it is convenient to use instead 
of (0) an equivalent set of equations for the density, and the local velocity of the condensate. The equations of motion 
are obtained by writing ip in terms of its amplitude / and phase (j), ip = /e**^- The number density is given by n = /^, 
while the velocity v is given by v = hVcjy/m. By inserting -0 = /e*"^ into (|]) and separating the equation into real 
and imaginary parts one obtains the two equations 



din f^. 



-V • if'Vcj,), (16) 



which is the equation of continuity, 



and 



dt m 



^+V-(nv) = 0, (17) 



We eliminate the phase variable by taking the gradient of dlq), using v — hWip/in. The resulting equation is written 
as 

m-^^^S/{Sfi + -mv^), (19) 

where 

5^l = V + Uon-- -V^V^-Mo- (20) 

Since it is the gradient of (5/i, which enters the acceleration equation (|l9|), we are free to subtract a constant from 
5 11. We have in (GO) chosen to subtract the value of the equilibrium chemical potential /ig, which implies that Sfi is 
zero in equilibrium, i. e. under stationary conditions. The equation Sfi — is the time-independent Gross-Pitaevskii 
equation. 

In the Thomas-Fermi approximation, the kinetic energy term is neglected, and the equilibrium density is therefore 
given by 

noUo + V{x,y,z)^ fiQ. (21) 

Within the Thomas-Fermi approximation one also neglects the kinetic energy term involving Sn in the expression ( |20| ) 
for (5/z. This yields 

Sfj. = Uo6n. (22) 

The Thomas-Fermi approximation for the modes should be a good one provided the wavelength of the mode is large 
compared with the healing distance S. With these approximations we may readily linearize the equations (O) and 
(|l9|), and eliminate Sfi by means of (|22|). The result is 

5ti 
m-g^ = C/oV • (noVSn). (23) 



If we only consider oscillations with the time dependence Sn ex exp(— zcji), the differential equation ( |23| ) simplifies to 

-u!'^6n=—(Wno-V6n + no^^6n). (24) 

m 

We investigate the surface modes in a two-dimensional configuration with the linear ramp potential considered in 
the previous section. In the y- and z-directions there is translational invariance, and therefore the solution must have 
the form of plane waves for these coordinates. We denote the wavenumber of the mode by q, and take the direction of 
propagation to be the z-axis. In the Thomas- Fermi approximation the condensate density in equilibrium, no, is then 
given by rig (a;) — —Fx/Uq for a; < 0, while it vanishes for x > 0. It follows that (p4[) has a solution of the form 



Sn^Ae"'""'', (25) 

which describes a wave propagating on the surface, and decaying exponentially in the interior. Since ( |25| ) satisfies 
V^(5n = 0, while the gradient of the equilibrium density is given by (— _F/[/oi 0,0), we obtain by inserting (E5^ into 
( P4| ) the dispersion relation 

w^ = —q. (26) 

ra 

This has the same form as for a gravity wave propagating on the surface of an incompressible ideal fluid in the presence 
of a gravitational field g — F/m. 



The solution (25) is however not the only one which decays exponentially in the interior. To investigate the solutions 



to ( |24| ) more generally we insert a function of the form 

(5n = /(ga;)e«"+'«^ (27) 



and obtain the following second order differential equation for f{y), 



y-^ + (2y + l)-^ + il-e)f^0, (28) 



where e = uP' /gq. By introducing the new variable z = — 2y one sees that Eq. (|8|) becomes the differential equation 
for the Laguerre polynomials L„(z), provided e — 1 = 2n. We have thus obtained the general dispersion relation for 
the surface modes 

p 
iJ = —q{l + 2n), n = 0,l,2--- (29) 

m 

with the associated density oscillations given by 

5n{x, z, t) = ALn{-2qx)e'i^+"^''^''^\ (30) 

with A being an arbitrary constant. 

To make contact with Stringari's calculation 0], we note that he found the dispersion relation of modes in an 
isotropic harmonic trap to be given by oj'^ = ujq[1{1 + 2n) + 3n + 2n^], where I is the angular momentum quantum 
number and n the radial one, which gives the number of nodes in the radial direction [pj. For I much greater than 1, 
the dispersion relation becomes cj^ = ^^oUl + 2n). The wavenumber of the mode at the surface of the cloud is given 
by 9 = l/R, and therefore the dispersion relation is w^ — cjgg_R(l + 2n), which is precisely the same as the result (p9) 
we obtained above, since the force due to the trap at the surface of the cloud is simply F = ^qR per unit mass. For 
large values of I it is thus a good approximation to replace the harmonic oscillator potential by the linear ramp, as one 
might expect since the surface modes are concentrated within a distance of order R/l from the surface. It should be 
noted that the n = mode frequencies for the plane surface with a linear ramp potential agree with the frequencies 
of the nodeless radial modes (corresponding to n = 0) for a harmonic trap at all values of I. For modes with radial 
nodes (n 7^ 0), the two results agree only when I is much greater than n. 

IV. SURFACE MODES AT SHORTER WAVELENGTHS 

When the healing length 6 is not negligible compared with the wavelength, there are corrections to the dispersion 
relation. In the case of gravity waves on the surface of a liquid, modes at shorter wavelengths are affected by the 
surface tension, and the dispersion relation is given by cu^ = gq + aq^ / p, where a is the surface tension, and p is 
the mass density of the fluid. We now explore modes on the surface of a Bose-Einstein condensed cloud at shorter 
wavelengths, and we shall show that there are contributions to w^ of order q'^ \n{l/q6) which may be understood in 
terms of the effective surface tension introduced in the discussion of static properties. 

The basic problem is to solve Eq. dl3) including the quantum pressure term in the expression for the chemical 
potential. Rather than attacking the problem directly, which leads to two coupled second-order differential equations, 
we shall adopt a variational approach, which will allow us to calculate the leading corrections to the Thomas-Fermi 
result for the mode frequencies for small q. 

In order to determine the dispersion relation of surface modes at shorter wavelengths we employ a trial wave function 
that allows us to calculate the total energy in terms of two variables which describe the displacement of the surface 
and the local velocity, respectively. In terms of these variables the energy functional assumes the form of that of a 
harmonic oscillator, from which we may extract the frequency as a function of q. 

The trial wave function is motivated by the solution found above in the Thomas-Fermi approximation. To lowest 
order we may describe the motion of the surface in a traveling wave by modifying the ground state wave function in 
two respects. First one shifts the spatial variable, thereby allowing for displacements of the surface, and one introduces 
a phase factor to take into account motion of the particles. Explicitly, the wave function is given by 

Tpix, z,t) = TpTpix — A(a;, z, t)) expicj) « il^Tpix) — ^{x, z, t)%j/rpp{x) + i(j)il)TF{x)- (31) 

Here A(a;, z, t) = ^0 exp qx cos{qz — tot), with ^0 constant, and 4> — {mvo/hq)e''^ si'ii{qz — tot), vq being the amplitude 
of the velocity at the surface of the cloud, given by 

V = = volsiniqz — cut), 0, cos(qz — ujt)). (32) 

m 

Let us now turn to the more general case. We expect that at frequencies small compared with the characteristic 
frequency ~ Ti/2m5'^ associated with adjustments of the density profile in the region within ^ J of the surface, the 



density profile will be able to adjust essentially instantaneously to its equilibrium form corresponding to the local 
number of particles per unit area, even if the Thomas-Fermi approximation is not valid. The real part of the wave 
function in the vicinity of the surface is thus of the equilibrium form, but with a possible translation perpendicular to 
the surface. We shall therefore use a trial function which has the same form as for the Thomas-Fermi case, but with 
the equilibrium Thomas-Fermi wave function replaced by the exact one. From this we shall calculate the energy, and 
evaluate oscillation frequencies. We write the wave function in the form 

V- = V'o + S^P, (33) 

where S'ip is the part due to the oscillation. For the present purposes it is simplest to consider a standing wave, and 
therefore we adopt the following form 

5^j{x, t) = (-e(t)Vo'(a;) + iM^)Mt)Vcos{qz). (34) 

The x-component of the velocity of the surface is given either in terms of the time derivative of the surface displacement, 
or in terms of the x-derivative of the phase of the wave function. This leads to the consistency condition ^ = Tiq(j)a/m. 
We use the trial function ( |34| ) to evaluate the energy functional 



E ^ dv 



'^iVVr + nr)|^(r,t)|2 + i[/oWr,i)|4 



(35) 



to second order in ^ and ^. The zeroth order term gives rise to an unimportant constant, while the first order term 
vanishes, since the trial function is a solution of the Gross-Pitaevskii equation. The interesting physics is contained 
in the second order term E^"^' . The ^^-part of the kinetic-energy contribution contains an integral over x of the form 

/>oo 

dx(V'o'e'?")'(Vo'e«")' = - / dxe^^-VolV-o'" + 2gVo" + ^zVo')- (36) 



The third-order derivative in (^6|) is eliminated by differentiating the Gross-Pitaevskii equation, 

- ^C = -V'i^o - VxPo' - 3t/oVoVo'. (37) 

The last two terms on the right hand side of ( pTf ) yield contributions that are canceled by those coming from the 
potential and interaction energies in (|35|). The remaining terms may be combined using V' — F and partial integration, 
and yield for the energy per unit area 



i-2 r /-oo i-2 i^ca 

i^(2) ^^ Fq dxe^^^^o^ + — dx (g^e^?-^^^) 

4 7-00 m J-oo 



""^^ ' dxe^^^Vo'. (38) 



This is of the same form as for an harmonic oscillator, E'^^^ — \Ci^^ + 'k^'iS,^ where Ci and C2 are constants, and 
the frequency is given by uj"^ — C1/C2. Integrating by parts the term involving _F, we obtain the final result 

C.2 = ^g-K^/(9), (39) 

where the dimensionless quantity I{q) is given by 

_ J^^dxe^i^i^o')^ 
q^ J_ dxe'^'i^'ipo 

The first term in (B3) gives the frequency of the surface mode in the Thomas-Fermi approximation, while the second 
term involving I(q) is a correction term. 

We may evaluate the leading long- wavelength corrections to the dispersion relation I{q) for q <^ 1/(5 by splitting 
up the range of integration into two regions 

/ — L /'oo 

dxe2«-(V'oT+ / dxe2«-(V^oT, (41) 

-00 J — L 



where 5 ^ L <C I/?. In the first of the two integrals we may use the Thomas-Fermi approximation -00 = \/^Fx/Uo, 
since L^ S. In the second we may replace exp2qx by 1, since ql <$^ 1. By using partial integration one finds 

/ dx{i;o')^ ^W rf^V'o'Vo (42) 

</ — J_j U i/ — J_j 

since the value of ipoip'o at a; = — i may be evaluated in the Thomas-Fermi approximation, resulting in —F/2Uq. The 
integral on the right hand side of (|4|) occurs in calculations of the kinetic energy associated with the surface region, 
and it has been evaluated previously. Its asymptotic form is (F/4{/o) ln(L/1.776^), which corresponds to the kinetic 
energy contribution in Thomas-Fermi theory cut off at a distance 1.776^, where the coefficient was determined by 
numerical integration. We obtain consequently 

^xe^-(W)^ = 4^(2-ln2,X-ln,) + ^ln^, (43) 



where 7 (~ 1.778) is the Euler constant. The denominator in (^0|) is evaluated in the Thomas- Fermi approximation, 
resulting in F/AUq. We finally obtain 

I{q) ~ - \iiq5 + ln(eV3.5527) ~ - InqJ + 0.15. (44) 

Thus we find 



u^ 



-q+^[-liiqS + 0.15]. (45) 



m 



The qualitative behavior is easy to understand by analogy with the surface tension contribution to the frequency of 
gravity waves on the surface of a fluid, where the contribution to w^ is aq^ / p. In the present problem, the effective 
surface tension depends logarithmically on the length scale, which is given by the wavelength. The q'*-dependence 
exhibited by ( [45[ ) then results from dividing g^ by the effective density in the region where the fluid is moving, this 
being of order the fluid density at the distance 1/q from the surface, or F/U^q. Thus one sees that the surface tension 
is very weakly dependent on the trap parameters and the atomic scattering length, which occur only in the logarithm. 
This result agrees to the order indicated with the result of using the full solution -00 {x) to the (equilibrium) Gross- 



Pitaevskii equation in evaluating the expression (40). In a recent paper. Fetter and Feder |g] analyzed the corrections 
to the excitation frequencies in the Thomas-Fermi limit 0] for an atomic cloud confined by a spherically symmetric 
trap. By employing the matching conditions of boundary-layer theory |7|, they were able to demonstrate that the 
leading correction to the Thomas- Fermi limit is of order i?""*, where R is the Thomas-Fermi radius, while terms of 
order R~'^ ln_R were found to be absent. Repeating their analysis for the different geometry which we are considering, 
we have explicitly verified that the matching conditions allow for the presence of terms of order q'^ Inq, which we have 
found in the present paper. 

This result for the surface mode frequency can be obtained in a more rigorous fashion by a variational approach. 
The resul ting equations of motion of ^ and 0o take the form of that of a classical harmonic oscillator with frequency 
given by (p9[). The details of this calculation are described in Appendix A. 

We thank Emil Lundh for helpful contributions. 

APPENDIX A: VARIATIONAL APPROACH 

The Gross-Pitaevskii equation may be derived from the variational principle 

6fdtL = 0, (Al) 

where 



Here E is the energy functional given by equation (|35[). For tp we adopt our ansatz, (|33|), and the equations of motion 
for ^ and 00 may then be determined from the variational principle. To second order in ^ and 0o the Lagrangian is 



L = -hi^<j>o-M)-iE^^^+E^'^), 



(A3) 



where E^^' is the ground-state energy obtained by inserting ip = ipo m equation (pa), while E^'^' is given by 



in terms of the integrals 



£;(2) = i! 
4 



h 



Fqh + —h 
m 



2„2 



h'q 
Am 



»0^2, 



dxe2«-(^2) ^ 



(A4) 



(A5) 



and 



dxe2«^ Vo- 



ls = / dxiq'e^^-^, ). 

>/ — OO 

The resulting equations of motion for ^ and <pQ are 

77 A '/'o - 7:aFqh + —h) = 
2 2 m 



and 



0^1?+ V^ ^200=0. 

2 2m 



(A6) 



(A7) 



(A8) 



(A9) 



By partial integration of ( [A5| ) we see that /i = —(7/2, and subsequently, by eliminating 00 between the last two 
equations, one obtains the result (|39|). 
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